Calculation of all elements of the Mueller matrix for scattering of light from a 
two-dimensional randomly rough metal surface 
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We calculate all the elements of the Mueller matrix, which contains all the polarization properties 
of light scattered from a two-dimensional randomly rough lossy metal surface. The calculations are 
carried out for arbitrary angles of incidence by the use of nonperturbative numerical solutions of 
the reduced Rayleigh equations for the scattering of p- and s-polarized light from a two-dimensional 
rough penetrable surface. The ability to model polarization effects in light scattering from surfaces 
enables better interpretation of experimental data and allows for the design of surfaces which possess 
useful polarization effects. 
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When light undergoes scattering from a surface, the 
scattered light carries a great deal of information about 
the statistical properties of the surface in its polarization. 
Even when the structures in question are sub- wavelength 
and beyond the imaging limit, polarized optical scatter- 
ing can be employed to detect and distinguish between 
material inhomogeneities, particles, or even buried de- 
fects and the roughness of both interfaces of thin films Q] . 
These techniques are already in use in the semiconduc- 
tor industry, and such techniques could become impor- 
tant for surface characterization of photovoltaic materi- 
als and nanomaterials [H [3]. Many biological materi- 
als are optically active, meaning that polarimetric mea- 
surements may be applied for characterization of biolog- 
ical or hybrid [3] materials and even for the search for 
extra-terrestrial life 5j . Surface patterning has also been 
proven as a method for creating optical components with 
interesting polarization properties [5] . 

To extract this information from experimental data, 
one has to be able to model polarization effects [7]. The 
ability to calculate the polarization of the scattered light 
also opens the door to the possibility of designing surfaces 
that produce specified polarization properties of the light 
scattered from them. 

All the information about the polarization proper- 
ties of light scattered from two-dimensional surfaces is 
contained in the Mueller matrix [MTU], Yet, very few 
calculations of the elements of this matrix for a two- 
dimensional randomly rough surface have to date been 
carried out by any computational approach, largely be- 
cause calculations of the scattering of light from such 
surfaces are still difficult to carry out pTMl4) . An excep- 
tion [15] is a calculation of the Mueller matrix for two- 
dimensional randomly rough perfectly conducting and 
metallic surfaces characterized by a surface profile func- 
tion that is a stationary, zero-mean, isotropic, Gaussian 
random process, defined by a Gaussian surface height 
autocorrelation function. These calculations were carried 



out by a ray-tracing approach on the assumption that the 
surface was illuminated at normal incidence. In this work 
it was also shown that due to the assumptions of normal 
incidence and the isotropy of the surface statistics, the 
elements of the corresponding Mueller matrix possess cer- 
tain symmetry properties. Subsequently Zhang and Ba- 
har [TB] carried out an approximate analytic calculation 
of the elements of the Mueller matrix for the scattering 
of light from two-dimensional randomly rough dielectric 
surfaces coated uniformly with a different dielectric ma- 
terial. 

In this Letter we report the first step toward realizing 
the possibilities mentioned in the opening paragraphs. 
We present an approach to calculating, for arbitrary an- 
gles of incidence, all the elements of the Mueller matrix 
for the scattering of light from a two-dimensional ran- 
domly rough metal surface. It is based on nonperturba- 
tive numerical solutions of the reduced Rayleigh equa- 
tions for the scattering of p- and s-polarized light from a 
two-dimensional rough penetrable surface [13l [17] . 

The system we study consists of vacuum in the re- 
gion X3 > C( x ||)i where Xy = (x\,X2,0), and a metal 
whose dielectric function is e(uj) in the region X3 < C( x ||)- 
The surface profile function C( x ||) i s assumed to be a 
single- valued function of xy that is differentiable with re- 
spect to xi and X2, and constitutes a stationary, zero- 
mean, isotropic, Gaussian random process defined by 
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( xy — x|| ^ . The angle brackets 



here and in all that follows denote an average over the en- 
semble of realizations of the surface profile function, and 
1/2 

S = (£ 2 (x||)) is the rms height of the surface. Each 
realization of the surface profile function was generated 
numerically by the filtering method used in Refs. [T4]I18]. 

We begin by writing the electric field in the vacuum 
region x% > C( x ||) as the sum of an incident and a scat- 
tered field, E(x,i) = [EW(x|w) + E( s )(x|cj)] exp(-iu;t), 
where 
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(^1,^2,0), the unit polarization vectors are 
ei (k||) = k|| x x 3 , 
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9||X 3 ], er(q||) = q|| xx 3 , 
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, with Reao((?||) > 0, 



(w/c) 2 - q 

Imao(9||) > 0. A caret over a vector indicates that it is 
a unit vector. In terms of the polar and azimuthal angles 
of incidence (9q,4>q) and scattering (# S ,0 S ), the vectors 
ky and qy are given by k|| = (ui/c) sin # (cos </>o, sin0 o , 0) 
and qn = (w/c) sin# s (cos0 s , sin</> s , 0). 

A linear relation exists between the amplitudes 
(qy) and £g (ky), which we write in the form (a = 
p,s, P=p,s) 



^(qil^E^Mqillki,)^ >(k,|). 

It was shown by Celli and his colleagues [T7j that the 
scattering amplitudes i? ai g(qy |ky) satisfy the matrix in- 
tegral equation (the reduced Rayleigh equation) 
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with R pp and i? ps forming the first row of the matrix R, 
where 

/ (7lQ||)=y d2x \\ exp[-i7C(x||)]exp(-iQy -xy), (3) 



and a(p\\) — e(u))(us/c) 
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with Rea(py) > 0, 



Ima(py) > 0. The matrices A/*± (p|| |q||) are given by 
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These equations were solved by the method described in 
detail in [13]. First, a realization of the surface profile 
function on a grid of N 2 points within a square region of 
the X\X2 plane of edge L. In evaluating the qy-integral in 
Eq. ([2]) the infinite limits of integration were replaced by 
finite ones: |qy | < Q/2, and the integral was carried out 
by a two-dimensional version of the extended midpoint 
rule |19j using a grid in the q\q2 plane that is determined 



by the Nyquist sampling theorem and the properties of 
the discrete Fourier transform. The function /(7IQ11) 
was evaluated by expanding the integrand in Eq. ([3| in 
powers of C( x j|) an d calculating the Fourier transform 
of C™(xy) by the Fast Fourier Transform. The resulting 
equations were solved by LU factorization. 

The scattering amplitudes R a p(q\\ |ky) play a central 
role in the calculation of the elements of the Mueller ma- 
trix. In terms of these amplitudes the elements of the 
Mueller matrix, M, are [20] 

Mn = C (|i? pp | 2 + |i? sp | 2 + |i? ps | 2 + |i? ss | 2 ) 
M 12 = C (|i? pp | 2 + |i? sp | 2 - |i? ps | 2 - |i? ss | 2 ) 
M13 = C (R pp Rp S + R sp Rts + Rp S R PP + R ss R* p ) 
M14 = iC (i? pp i?p S + R sp R* s — R ps Rp p — i? ss i?* p ) 
M 2 i = C (|i? pp | 2 - |i? sp | 2 + |i? ps | 2 - |i? ss | 2 ) 
M22 = C (|i? pp | 2 - |i? sp | 2 - |i? ps | 2 + |i? ss | 2 ) 
M23 = C (i?p P i?* s — R SV R* SS + Rp S R pp — R ss R* p ) 
M24 = iC (i?p P i?* s — R sp R* s — Rp S Rpp + i? ss i?* p ) 
M 3 i = C (RppR^p + RgpRpp + RpsR^ s + RssR^s) 
M32 = C [RppR*p + RspRpp — RpsRs S — RssRps) 
M33 = C (i? pp i?* s + RspRp S + RpsR^p + RssRpp) 
M34 = iC (i? pp i?* s + R sp Rp S — Rp S R* p — R ss Rpp) 
M41 = — iC (i?ppi?* p — R S pRpp + Rp S Rl s — RssRps) 
M42 = — iC {RppRlp ~ RspRpp — RpsRts + RssRps) 
M43 = — iC (i? pp i?* s — RspRps — RpsRfp — RssRpp) 
M44 = C (i? P pi? ss — RspRp S — Rp S R S p + RssRpp) 

where 

(271-0 ) cosfto 

and L 2 is the area of the plane x 3 = covered by the 
rough surface. 

As we are concerned with scattering from a randomly 
rough surface, it is the average, (M) , of the Mueller ma- 
trix over the ensemble of realizations of the surface pro- 
file function that we seek. In evaluating an average of 
the form {R a pR* s ) we can write R a p as the sum of its 
mean value and its fluctuation about the mean, R a p = 
{Rap) + (R a /3 — (Rafs)) ■ We then obtain the result 

(R a/ 3R* 7S ) = (R a p)(R*s) + {{RapR^t) - {Rafs){Rj5))- 
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The first term on the right hand side of this equation 
arises in the contribution to an element of the ensemble 
averaged Mueller matrix from the light scattered coher- 
ently (specularly) ; the second term arises in the contribu- 
tion to that ensemble averaged matrix element from the 
light scattered incoherently (diffusely). It is the latter 
contribution, (M) lncoh , that we calculate. 

We have calculated in this way the 16 elements of the 
Mueller matrix when light of wavelength A = 457.9 nm 
is incident on a two-dimensional randomly rough silver 
surface whose dielectric function at this wavelength is 
e(lu) — —7.5 + i0.24 J2TJ. The roughness of the surface 
is defined by a surface height autocorrelation function 
W(|x|||) = exp(— Xy/a 2 ), where a = A/4 and the rms 
height 5 — A/40. For the numerical parameters we used 
L = 25A and N x — 319 which implies that Q = 6.4(w/c). 
For these parameters, and when the metal is assumed to 
be non-absorbing [Ime(w) = 0], our simulation approach 
conserved energy within a margin of 1% or better. More- 
over, the calculated Mueller matrices were found to be 
physically realizable and therefore self-consistent by the 
method of Ref. [22] • 

The results presented in Fig. [I] were obtained for an- 
gles of incidence (#o:0o) = (2°, 45°), i.e. for (essentially) 
normal incidence. The first thing to notice from Fig. [T] 
is that the individual matrix elements possess the sym- 
metry properties predicted by Bruce 15J. The elements 
of the first and last column are circularly symmetric; 
each element of the second and third columns is invari- 
ant under a combined 90° rotation about the origin and 
a change of sign; and the elements of the second column 
are 45° rotations of the elements of the third column in 
the same row [25] . Note that the elements (M 31 ) incoh , 
( M 40 incoh , (M u ) incoh , and ( M 24> incoh are zero to the 
precision used in this calculation. However, simulations 
indicate that this does not hold for anisotropic surfaces. 

The results presented in Fig. [3] were obtained for angles 
of incidence (#o,</>o) = (25°, 45°), and display some in- 
teresting features. The elements (Mn) incoh , (^22) incoh , 
and (M 3 z) incoh contain a (weak) enhanced backscattering 
peak at qy = — ky (Fig. [2]). The element (^44) incoh ap- 
pears to have a dip in the retroreflection direction. This 
dip is not present in the results of a calculation based on 
small-amplitude perturbation theory to the lowest (sec- 
ond) order in the surface profile function, and is there- 
fore a multiple scattering effect, just as the enhanced 
backscattering peak is. In contrast to what was the 
case for normal incidence, the elements (Msi) i . and 
(M 24 ) incoh are no longer zero. 

If we denote the ensemble average of the contribu- 
tion to a normalized element of the Mueller matrix 
from the light that has been scattered incoherently by 
rriij = (Mij) incoh / (Mn) incoh , we can estimate the or- 
der of magnitude of the Mueller matrix elements by 
calculating the quantities sy = (|mij(qy)|) , where 



(/(q||)) q = / d2 <?|| f(<i\\)/^(^/c) 2 , and the integral over 
q|| is taken over the circular region < q\\ < uj/c. It was 
found that s u , s 22 , s 23 , s 32 , s 33 , s 44 are of 0(1); si 2 , 
S13, s 2 i, s 34 , s 42 , s 43 are of O(0. 1); and s u , s 24 , s 31 , s 41 
are of 0(0.01). These results are only weakly dependent 
on the polar angle of incidence 9 , for the values of 8 
assumed in this study. 

In conclusion, we have presented a new approach to 
the calculation of all sixteen elements of the Mueller ma- 
trix for light scattered from a two-dimensional, randomly 
rough, lossy metal surface, for arbitrary values of the po- 
lar and azimuthal angles of incidence. It is based on 
a rigorous numerical solution of the reduced Rayleigh 
equation for the scattering of p- and s-polarized light 
from a two-dimensional rough surface of a penetrable 
medium, that captures multiple-scattering processes of 
all orders. The results display multiple scattering effects 
in certain matrix elements, such as an enhanced backscat- 
tering peak in the retroreflection direction, or an unex- 
pected dip in the same direction. The matrix elements 
also display symmetry properties that, for normal inci- 
dence, agree with those predicted by Bruce [15] . 

The physical implications of the approach and results 
of this Letter point to better understanding of the po- 
larimetric properties of random surfaces. Such knowl- 
edge may be critical for improved photovoltaic and re- 
mote sensing applications and has the potential of being 
used to engineer surface structures with well-defined po- 
larization properties of the light interacting with them. 
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FIG. 1. (Color online) Color-level plots of the contribution to the Mueller matrix elements from the light scattered incoherently 
as functions of qi and qi for angles of incidence (8o,(j>o) = (2°, 45°). An ensemble consisting of N p = 10 000 surface realizations 
was used in obtaining these results. The elements, {Mij) incoh = 1,2,3,4), are organized as a matrix with (Mii) incoh in the 
top left corner; (Afi2) incoh top row and second column, etc. The white spots indicate the specular direction in reflection. 
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f ps in Bruce's Eq. (1) corresponds to our R ps and R sp , 
respectively. 
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FIG. 2. The incoherent contribution to the diagonal Mueller 
matrix elements in the plane of incidence (parameters as in 
Fig. [3j|. The vertical dotted line indicates the backscattering 
direction. 
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FIG. 3. (Color online) Same as Fig. [I] but now for angles of incidence (9o,cj>o) = (25°, 45°). 



